function S = assemble_Stiff(dofs,V,T,nt,d,MassInt,Stiff_UU,Stiff_UV,Stiff_VV,Indx1_cell,Indx2_cell)
% prepare spaces for nonzero entry
m = 0;
for k=1:nt
    md = (d(k)+1)*(d(k)+2)/2;
    m = m + md^2;
end
Indx1 = zeros(m,1);
Indx2 = zeros(m,1);
S_stiff = zeros(m,1);
% begin to start
pos = 1;
  for tri = 1:nt
      degree = d(tri);
      v1=V(T(tri,1),:);v2=V(T(tri,2),:);v3=V(T(tri,3),:);
      x13 = v1(1) - v3(1);  %a
      x23 = v2(1) - v3(1);  %-c
      y23 = v2(2) - v3(2);  %d
      y13 = v1(2) - v3(2);  %-b
      J = x13*y23 - x23*y13;  
      a1 = (x23^2 + y23^2)/J^2;  
      a12 = (x23*x13 + y23*y13)/J^2; 
      a2 = (x13^2 + y13^2)/J^2; 
      % get the dof map
      loc_dof = get_tri_dof(dofs,tri);      
      L = length(Indx1_cell{degree});
      loc_pos = pos:pos+L-1;
      Indx1(loc_pos) = loc_dof(Indx1_cell{degree});
      Indx2(loc_pos) = loc_dof(Indx2_cell{degree});
      S_stiff(loc_pos) = (a1*Stiff_UU{degree} - a12*Stiff_UV{degree} + a2*Stiff_VV{degree})*J;
      pos = pos + L;
  end
dim = max(max(dofs));
S = sparse(Indx1(1:(pos-1)),Indx2(1:(pos-1)),S_stiff(1:(pos-1)),dim,dim);